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Abstract 

We perform full-scale numerical simulation of instability of weakly nonlinear 
waves on the surface of deep fluid. We show that the instability development 
leads to chaotization and formation of wave turbulence. 

We study instability both of propagating and standing waves. We studied 
separately pure capillary wave unstable due to three-wave interactions and pure 
gravity waves unstable due to four- wave interactions. The theoretical description 
of instabilities in all cases is included into the article. The numerical algorithm 
used in these and many other previous simulations performed by authors is de- 
scribed in details. 
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1. Introduction. 



Stationary propagating waves on the surface of deep heavy ideal fluid are 
known since the middle of nineteenth century. Stokes (see, for instance Jl]]) in 
1847 found the solution of the Euler equation in the form of trigonometric series. 
For the shape of surface rj(x, t) he obtained 



1 3 , 
co^{kx - ojt) + -jx cos{2(fcx - ii)t)\ + -jjL cos{3(fcx - iot)] + . . . 

2 8 



• (1) 



?7(jc, t) = a 

Here we introduced steepness n and frequency co 

(2) 



fx = ka, co= V^|l + + + ■ • • j 



Stokes found two algorithms for calculation of all terms in series (HI) and © 
See Sretenskii f^). Convergence of these series has been proven by Nekrasov [J, 
in 1921. Another proof was found by Levi-Chivita [5J. 
It is known since 1965 ^ that stationary waves on surface of deep water 
are unstable. The theory of instability [[Zl Isl Isl UM was developed for waves of 



small amplitude in the limit ^ 0. A history of this question is described in the 
article [11]. In the present paper we study instability of stationary waves numer- 
ically by direct solution of the Euler equation describing a potential flow of ideal 
fluid with free surface. This approach has two following important advantages. 
Firstly, by numerical simulation we can study waves with finite amplitude. This 
is mostly a program for future, in this paper we shall study only waves of small 
amplitude. The second advantage is not less serious. By the use of numerical 
simulation we can study not only linear, but also nonlinear stage of instability 
development. Even in integrable systems like the NLSE, analytical study of the 
monochromatic wave is a very nonlinear problem which could be solved only by 
methods of algebraic geometry. In more realistic models development of a non- 
linear theory of modulation waves instability is a hopeless problem. In a long 
run we have to expect that the instability will lead to formation of a stochastic 
wave field described by a kinetic equation for squared wave amplitudes and for- 
mation of Kolmogorov-Zakharov (KZ) spectra, governed by the energy flux in 



high wave numbers [Il2n 



The article is organized as follows. Chapters 2 and 3 are devoted to analytical 
theory of stability of weakly nonlinear stationary waves. To develop this theory 
we use Hamiltonian formalism as the most compact and suitable. We start with 
presenting the Euler equation of ideal fluid with free surface in the Hamiltonian 
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form. Surface tension is included in the Hamiltonian. In the presence of surface 
tension the dispersion relation is 



cok = yjgk + crk^, 

here cr - surface tension coefficient. 

Wave vectors of small-amplitude stationary waves are solutions of equation 

ojk = ck. (3) 

This equation has two solutions 



K2 = ^ —, (4) 

if c > Co, where Cq = For water cq - I2cm/sec. In a generic case 

c Co stationary waves comprise a complicated four-parameter family. However, 
in the limiting case c » cq one can split it in two periodic families of "pure 
gravitational" and "pure capillary" waves. 

The Stokes wave is "pure gravitational". Now ki = g/c^ and capillary ef- 
fects could be neglected. In the "pure capillary" case ^2 = c^/cr and effects of 
gravity can be neglected. All stationary waves on the surface of deep fluid are 
unstable. However, instability of short capillary waves and long gravity waves 
have different flavor. They described by different "efficient Hamiltonians". The 
case of "pure capillary" waves is the simplest. The instability can be studied if 
framework of the Hamiltonian contains only quadratic and cubic terms. This is 
a subject of Chapter 2. A situation is more complicated for gravitational waves. 
In this case terms of fourth order must be included in the Hamiltonian. Then 
one have to exclude the cubic terms by a proper conformal transformation. As 
a result we get so-called "Zakharov equation" [lOJ. In the framework of this 
equation the problem of the Stokes wave stability can be solved exactly. This is 
a subject of Chapter 3. 

In Chapter 4 we give a detailed description of the numerical code which we 
used for solution of the Hamiltonian Euler equation. This code was used in 
many papers but never was described in details [13l [lH [isl [l6l,[l7, 18. 19. 201. 



We should stress that in our numerical experiments we worked with the Euler 
equation written in "natural variables". These equation are not as good for direct 
analytical study as they good for implementation of numerical method. The 
structure of nonlinear parts of the Hamiltonian in "natural variables" is relatively 
simple and using the standard Fast Fourier Transform (EFT) method is quite 
possible. 
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In Chapter 5 we present our results on modeling of capillary wave instability. 
We show that an initial stage of instability is described pretty well by the linear 
analytical theory. Further development of instability demonstrates appearing of 
"secondary instabilities" and a tendency to formation of a chaotic wave field 
which should be described by statistical methods. 

In Chapter 6 we study instability of the Stokes wave. We show that this 
instability is mostly "modulational". In other words the wave remains quasi- 
monochromatic for a long time after development of the instability. 

Finally, in Chapter 7 we present first results on development of the standing 
waves instability. We show that this instability leads to fast isotropization of 
the wave field. This mechanism can be used in experiments for generation of 
isotropic wave field. 

2. Theory of decay instability 

In this section we develop the simplest version of the theory of stationary 
waves instability. This simple theory is applicable if triple- wave nonlinear pro- 
cesses governed by the resonant conditions 

o)k = cok^ + cok,, (5) 



are permitted. But before we briefly describe how the theory of surface waves 
can be embedded into the general Hamiltonian theory of nonlinear waves. 

Suppose that ideal incompressible fluid fills the space -oo < z < r}{f, t), 
here r = {x,y) — two dimensional vector. A flow is potential v = VO, hence 
hydrodynamical potential O satisfies the Laplace equation 



Let us define i/r = 01^-=,, and impose a natural boundary condition Phi, 
at z ^ -oo. It is known [7] that rj{f, t) and i/^(r, t) are canonically conjugated 
variables satisfying evolutionary equations 



dT] _6H dip _ 6H 
Here H = T + U — total energy of the fluid, consisting of kinetic energy 



k = k\ + k2. 



(6) 



AO = 0. 



(7) 




(9) 



— oo 
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and potential energy 

[/ = I ri^d^f+(T J ( Vl + (V77)2 - l)d2^. 

The Hamiltonian H in terms of 77 and iff is given by the infinite series 
H = H0 + H1+H2 + ... 

Here 

here ^i/r = V-Ai/r, 

1 r 

Thereafter we will neglect last term in (fT4] ). 

One can perform the symmetric Fourier transform 



^0 = ^ J + ^^7' + TiVr]f} dV, 
= ^ J 77{|V«AP-(W}dV, 



This is the canonical transformation. Equations dS) now take the form 
dT] _6H dil/ _ 5H 



Now 



Ho 
Hi 



Here 



L,-,- = (kik2) + MM, 



+\h + k2\ + \k2 + h\)-\ki\-\k2\ 



Equations ([8]) written for the Hamiltonian (11311141) read 



f] = kif/- (V(77Vi/r)) - k[T]kil/] + k{T]k[T]ki//])+ 
-1V2[,7%] + i 



+^V2[77V] + ^^['y'VV], (19) 



(/r = a-V^T] -gri-\ [(ViA)^ - - [kil^M-nH^ - [7lkil/]V^if/. 

These equations were considered for the first time in [12 ill . Equations (fT9l) are 
basic in our numerical simulations. To develop analytical theory of stationary 
waves instability we use equation (fT6l) . 

Let us introduce complex normal variables 



As far as ri_^ = rjl, if/ j^ = if/*^ (because these are Fourier transforms of real 
functions) we have 



= \i — + 'f'i = A/ - (21) 



In terms of equations (|T6l) turns to one equation 

% SH 

k 

Now 

//o = J a;^|a/d2^. (23) 

Then 

H, = Hf'U H'^'^\ (24) 

Here 

^(0,3) = 1 1 ^(O'l^ (a,-aj a^r + a*Ml al )5{t + + 4)d^d^id^2 (25) 

//f ''^ = ^ / ^Sa^''^^''^^'''^^^ ^ )^(^ - ^1 - 4)dM^id4. (26) 

In a similar way 

H2=Hf^'^+Hi'^'UHf^^\ ill) 



Only the last term in H2 is important for us 



^(2,2) = i r y(2.2) ^ ^* ^* + ^ ^ ^ 4)dM^id4d4. (28) 



Explicit expressions for Vl^'^i , Vli'^i , and vS'^i ^ are presented inlAppendix B 

kk\k2 kk{k2 kk[k2k^ 



Now everything depends on shape of function <jj{k). If resonant conditions ([5]) 



have real solutions, one can neglect H2 and even . Now equations takes 
a simple form 



4/' 



i<:i<:2 ' - kikk2 K2 



Thereafter we assume 

V^}'}\ = 0. (30) 

For surface waves this condition is satisfied. 

Equation (l29l) has a solution which can be treated as a stationary wave 

00 

= [an^~"'^'5{t - nto) + bne'"^'6(,l^ + n^)] • (31) 

n=l 

-> 

We put ai = s. Here /:o — an arbitrary wave vector. Coefficients a„ are presented 
by power series 

while coefficients look as follows 

K = s"^\bf + eU;^^ + ...). 

The frequency of the stationary wave is presented by a series in even powers of 
s 

^ = coiko) + e^Ai + £^A2 + ... (32) 

Now let us suppose s ^ 0. From we see that the first nonlinear correction 
to frequency is proportional to s^. Now we will show that the solution (|3T]) is 
unstable and the growth rate of instability is proportional to s. It means that all 
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nonlinear corrections to a„, b„, and A„ can be neglected, and one should look for 
a solution in the following form 

= £e-'"^'^°)'5(^ - ko) + ait)e-""-''^'6{k - + me-^^-'^^'Sik - 1I2), (33) 

where k\ + K2 = h). Then we linearize the equation (1291 ) and find that a and yS 
obey to the system of ordinary differential equations 

a = ie-'^sV/3\ V = V^^' , 

koi<lK2 

= ie~'^sVa\ A = co{ko) - to{Ki) - oj{k2) 
A general solution of equation (|34l) is 

a = aoo^-'^'^^^^ fi=/3oe^''^'^^^^'. (35) 

Here 



koK^i (34) 



7 = y e'lVp - (36) 
and ao, /3q are connected by relation 

~±yyo = isVa;. (37) 

Instability takes place if wave vectors /ci , K2 are posed near the surface (or the 
curve) 

oj{ko) = a)(Ki) + a)(K2), ^^^^ 

ko = Ki+ K2. 

The maximum of the growth rate y = £\V\ is reached exactly on the resonant 
surface, which is represented on Figure [1] The plot of growth rate (1361 ) on a 
discrete grid of wave numbers for In x In periodic box is given in Figure |2l 

3. Four-waves instability. 

The theory developed in the previous chapter is applicable to study of capil- 
lary waves instability. In the case of gravity waves resonant conditions (|38l) have 
no real solutions. It means that cubic terms in the Hamiltonian can be excluded 



by a proper canonical transformation given by power series [1221 12311 



a, = bf + bi'' + b? + ... 

(39) 

I- 



^ k k ' li 
b'? = a 



k 
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Figure 2: Growth rate for decay instability of monochromatic capillary wave with feo = (68, 0) 
and average steepness yu = 0.05 on the discrete grid. Upper panel: isometric projection; lower 
panel: contour of the surface at the level 10"^^. 



10 



Next terms b''P, Z??^ are presented inlAppendix B 

k k ' 



After canonical transformation the Hamiltonian is reduced to the form 

H = Ho + 62, (40) 
^2 = \J Tj;j^j^jblb*.b^bj^^5(k + h-k2- 4)dM^id^2d^3. (41) 



Explicit and complicated expression for T^^^^^^^j is given in 
is a homogeneous function of order three. 



Appendix B[ 



Now canonical variable Z?^ obeys so called "Zakharov's equation" IllOll 



Tii.uJ'lhH^^^ + ki-k2- k,)dkidk2dk3. (42) 



Equation (l42l) has exact solution 
bj^ = A6{k - 

ao = oj{h) + ]^T-,^\A\\ 



(43) 



-ko 



— - = —k^ 

kokokoko ^' 



This solution is nothing but the stationary Stokes wave. It gives right expression 
for two first terms in series (H) and (|2). 

Equation (1421) has also reach set of approximate quasi-periodic solutions. Let 
^1 and K2 are two arbitrary wave vectors. In the limit of small |Aip, \A2\^ equation 
has the following solution: 



bj- = Ai5{k- /?i)e"'"'' + A26ik - y?2)e 



-lQ2t 



+ ... 



(44) 



2 ^ (45) 

^2 = C0{K2) + r,^,,^JAi|2 + -Tg,\A2\\ 

Here r^,/2 = 7'/?i/?2,'?i'?2- Equations (l44|) . (|431) are valid if nonlinear terms in (05) 
are much less that linear. 

In the particular case K2 = -i<i, IA2I = IAi| solution (l44l) is just a standing 
wave. 

Both propagating wave (1431) and standing wave (|44l) are unstable. To study 
instability of propagating wave (|43]) we will look for a solution in the following 
form 

bj; = Ao6{k - 4)e"'"°' + aS{k - ko - ^e"'"" + /36{k - 4 + ^e-'"^'. (46) 
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Here 

Qi = oj{ko + i<) + 2T{kQ, 4 + k)\A\^, 



2 (47) 



By plugging (l44l) into (l42l) and linearizing over a and /3 we set system of ordinary 
differential equations, similar to (|34] ) 

$ = ^c-'^\Ao\^Ta*, (48) 
A = 2oj(ko) - co(ko - - co(ko + + 2(r^^ - T^^^^^ - r^„,^,)|AoP. 

Solutions of equation (l48l) are given by formulae (|35T ) where 



7= y|Anrp-^A2. (49) 

One can see that instability occurs in a vicinity of the curve A = 0. For waves on 
a deep water we can put k^ = i, k = xi + y j. Then the condition A = is reduced 
to the famous Phillips curve n2m 

[(1 + x)2 + [(1 - xf+y^f' = 2. (50) 

Here -5/4 < x < 5/4. The Phillips curve is plotted in Figure |3] The coupling 
coefficient T{ko, x), evaluated on the Phillips curve, can be presented in the form 

T(2kQ, ko + K,ko - it) = klf(x), x = 

/Co 

f{-x) = f(x) symmetric function. It is important to mention that /(5/4) = 0. 
This fact was first discovered by Dyachenko and Zakharov in 1994 [12511 . De- 
creasing of f{x) with growth of x means that the four-wave instability is mostly 
modulational, because the most unstable modes are concentrated at /c ^ 0. In 
this region 



Instability is concentrated inside the angle where 



^12 



KaKp < 0. 
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Figure 4: Growth rate for Phillips instability of monochromatic gravity wave with - (30, 0) 
and average steepness /i = 0.1 on the discrete grid. Upper panel: isometric projection; lower 
panel: contour of the surface at the level 10"^^. 
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The plot of growth rate (|49] ) on a discrete grid of wave numbers for In x In 
periodic box is given in Figure |4] Four-wave instability of propag ating Stokes 



waves was studied in details in many papers (e. g. [6, viliilS, 



SliilS). 



Let us study instability of by-harmonic wave (|46l) . We concentrate only on 
the case of standing wave K2 = -ki = k^. The standing wave is unstable due to 
a different mechanism. Firstly, each propagating wave composing the standing 
wave endure its own modulational instability. Then another instability appears. 
In this new type of instability we have simultaneous excitation of two waves with 
wave numbers +/?such that \k\ = |A;o|. The maximal growth-rate of this instability 

Here 

f = T^^^-u,,t-n = klf{cose), (52) 
-> 

and 6 is an angle between and k. 

4. Numerical Simulation Scheme. 

The problem of numerical integration of system of equations (fT9l ) is rather 
hard. One of the most important questions is what time integration algorithm to 
choose. According to the property of the equations it would be natural to develop 
a numerical integration scheme conserving the Hamiltonian. Let us follow the 



article fl27ll . One can introduce a discrete variation of Hamiltonian (fTTl) on one 
step on time H" — > H"*^ 

AH = H"^^ - H". (53) 

The Hamiltonian is a function of canonical variables 77 and iff. The discrete vari- 
ations of these functions on a time step are equal to 

A77 = 77"+! - 77", Aiff = iff"^^ - iff". (54) 



One can expand discrete variation AH via Arj and Ai{/ (that is done in Appendix A ) 



AH = H^Ail/ + H^Ari. (55) 

It is easy to see that Hr^ and are discrete analogues to the continuous variations 
f and f . 

or) 01// 
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One can demand conservation of Hamiltonian AiZ/r = during time step r. 
Obviously, this equality can take place if the following conditions are valid 



- n^. 



(56) 



In some sense this is a discrete analogue of Hamiltonian equations (fT6l) . Thus, if 
Hamiltonian variation (1531) is expanded via variations Ar] and At/^, it is possible 
to get equations (l56l) . 

As it was mentioned above it is more convenient to rewrite the equations in 
terms of Fourier harmonics. Using the results obtained in Appendix A (|A.2HA.7l ) 
one can obtain an implicit difference scheme 



n+l 
'k 



'k 



k k 



-\F (V, (?7"+i + 77")V(iA«+i + i/r")) - 
+\\k\F [(?7"+i + ri")k(r]"'-^kif/"^^ + T]"kif/")] 
+i\li\F [((.r]"^^f + (?7")^)V2((/."+i + i/r")] . 

2 \k\ V 'k 'ii 



+ 



Here F is the Fourier transform operator. 

It is useful to resolve linear part of scheme 



(57) 



(58) 



with respect to 77"^' and 
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if/"'^ . Let us denote nonlinear terms in right hand sides of these equations as: 



R 



-\F (V, (77"+i + T]")V(if/"^^ + i/^")) - 
-l\k\F [(7]"^^ + 77")%"+' + i/r")) + 
+\\k\F [(?7«+i + r]")k{T]"+^kif/"+^ + ri"ki/f")] - 
-l\k\^F [((77"+^)2 + {r]''f)kiil/"^^ + + 

= -|Ff|V.A""'r + |V«A"l') + 

|F + tf/")k(T]"^^ki//"^^ + r]"ktf/")] - 



+ 



1 ^ 

4f 



Using these notations discrete scheme can be written as follows 



r fi+l 



lA'r ' = E(k, t)t]'1 + A{k, T)if/Z + F(k, T)Rf^ + C{k, t)R 



5J1+1 



(59) 



(60) 



Here 



A{k,T) 
C{k,T) 

E{k,T) 



1 - >v 

1 + 



B(k,T) = 



Tk 



1 + 



D(^,T) = -TB{k,T), 



-^C{k,T), Fik,T) = -TEik,T). 

k 2 



(61) 



Thus, we get implicit (terms R"'^^ and RV'^ contain 77""^' and iff'l'^^) difference 

I Y k k 

scheme. The important feature of this scheme is that conservation of Hamilto- 
nian (fT3])- (fT4l) is embedded in it. 

The implicit numerical scheme (1601 ) can be solved by method of simple it- 



erations. Let us write this procedure for 77"^^'* and il/', 

k k 



n+\,s 



here s is an iteration 



number. Corresponding to (|60d one can get 



'k 



'k 



= A{k, T)7]'l + B{k, tWI + C{k, r)i?r''' + ^(^' ^)^r'''' 

- E{k, T)rfl + A{k, tWI + F{K T)i?r''' + C{K t)R;^'''; 



(62) 



Iterations continues until the desired accuracy of Hamiltonian conservation e 

is achieved. In most cases it is enough to follow the convergence of the relative 
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error 




n+l,s+\ 



,'1+1, * 



< 6. 



(63) 



When studying gravity waves this condition is equivalent to calculation of po- 
tential energy with desired accuracy. For the weakly nonlinear regime quadratic 
part of Hamiltonian is dominant, so the physical meaning of this condition is 
quite clear. 

This numerical scheme can be used for simulation of freely decaying waves. 
For simulation of turbulence we should introduce pumping and damping terms in 
the equations. It is possible to do that by different ways. We have applied split- 
step method, which is highly used in numerical simulation of pulse propagation 
in optic fibers. Let us consider a simple example. 

Let us suppose a linear damping with rate yi, in our model 



It is possible to take into account this damping without significant changes in 
calculations scheme. First one can obtain solution of equations (fT9l) without 
damping using iteration scheme described above. Let us denote this solution by 
i/r^+i . Second, the solution of the whole system of equations can be calculated by 
the next step 



It is worth to say that for weak turbulence simulation the most interesting part of 
spectrum is in the "inertial interval" where there are no damping or pumping at 
all. Even more, the nature of damping and pumping is not important. In this case 
influence of non-conservative terms can be described by such a rough scheme. 
As a bonus we have eliminated the restrictions on time step 

maxdr^Dr < 1, (66) 

unavoidable in the case of integration by standard Runge-Kutta methods. 
Pumping can be considered in a similar way. 



if/j^ = R.H.S.-yki/fj^. 



(64) 



= {jjf' exp(-y,r). 



(65) 
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5. Capillary waves 



In this section we briefly review our previous results published in [|l3ll and 
report new observations. System of equations ([T9) was simulated in the domain 
Lv = Ly = In. Surface tension coefficient cr = 1. Number of grid points 
was 512 X 512. A monochromatic wave of amplitude |a^J = 2 x 10"^, which 
corresponds to average steepness jj. = 0.05, was taken as initial conditions. Its 
wave number vector = (68, 0). All other harmonics were of amplitude |a^| ~ 
10~'^ and with random phase (Figure |5]). As was mentioned above resonant 




Figure 5: Decay of the monochromatic capillary wave. Initial conditions. Time f = 0. 

curve almost never passes through grid points (there are two non trivial points 
k = (0;0) and initial wave k = (ko,Oy, this process has a zero growth rate). 
Detailed picture of resonant curve on the grid in region with highest grid points 
density in the vicinity of the curve is shown in Figure |6l One can see, that some 
points are closer to resonant curve than others. 

In the beginning one can see growth of several harmonics as it is predicted 
in (|35l) and (l36l) . Different stages of the decay process are represented in Fig- 
ures iTlfTOl Time is given in periods of initial wave Tq. We represent isometric 
projection of the |a^p-surface and contour of this surface at the level 10"^^ (order 
of magnitude higher than background noise). The full picture of the fc-plane in 
the final moment of simulations is represented in Figure [H] Zoom of the most 
interesting region of the Figure [TT] (initial decay region) is represented in Fig- 
ure [12] One can see that although the amplitudes of the waves are stochastic, the 
spectrum is still strongly anisotropic. 

19 




Figure 6: Part of resonant curve for decay of monochromatic capillary wave with ko - (68,0). 
The different mismatches for different grid knots are clearly seen. 
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Figure 12: Decay of the monochromatic capillary wave. Zoom of the initial decay region. Time 
t = I44488r0. 
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6. Gravity waves. 



In the case of gravity waves on the surface of deep fluid the dispersion is the 
following 

= V^, (67) 

here g is a gravity acceleration. Here and further let us suppose ^ = 1 . 

In this case dispersion is of non-decay type conditions (I^Sl) have no real non- 
trivial solutions, and main process is four-wave scattering. Therefore one can 
make a substitution to eliminate third order terms corresponding to the decay 
process. This is the reason why we have to use Hamiltonian expansion up to 
forth order in the case of gravity waves. 

Let us consider the same initial conditions as in the case of monochromatic 
capillary wave decay, i.e. one monochromatic wave and random phase noise 
of small amplitude. The main processes involve large amplitude of initial wave 
most times. In this case one wave to three and invert processes are much weaker 
than scattering two waves with the same amplitude and the same wave vector to 
two other waves. 

Resonance conditions for such process are as follows 

cou^ + u)k2 = ^ojko^ ^1 -I- 4 = 24- (68) 

The resonant curve for this conditions is shown in Figure |3 

System of equations (fT9l) was simulated in domain = Ly = In, gravity 
acceleration ^ = 1. Grid size was equal to 512 x 512 points. As an initial con- 
ditions monochromatic wave of amplitude |a^J = 1.3 x 10"^ with wave number 

vector ^0 = (30, 0) was used. All other harmonics were of amplitude \a^\ ~ 10"'^ 
and random phase (Figure [13]). In the beginning one can observe exponential 
growth of several harmonics in the vicinity of resonant curve (detailed picture 
of resonant curve in the surroundings of the initial wave is shown in Figure [T?l). 
This is shown on Figure [T5l It is clearly seen that wave with wavevector (33, 2) 
has smallest mismatch and, due to that, grows. As we already know four-wave 
scattering growth rate has highest values in the vicinity of ho = (^0,0). Due to 
this the initial growth is lumped about the cross of resonant manifold. Different 
stages are represented in Figures [T6]|2n We represent isometric projection of the 
|a^p-surface and contour of this surface at the level 10"^^ (order of magnitude 
higher than background noise). The full picture of the fc-plane in the final 
moment of simulations is represented in Figure |22l Zoom of the most interest- 
ing region of the Figure |22] (initial instability region) is represented in Figure |23] 
One can observe still weak but observable downshift of the spectrum. 
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Figure 14: Gravity waves. Part of resonant curve. Different mismatch for different grid points is 
well seen. 
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Figure 15: Gravity waves. Growth of harmonics amplitude as a function on time. One can see, 
that harmonic k = (33, 2) is in almost exact resonance and the others are not. 




Figure 16: Instability of the monochromatic gravity wave. Growth of the harmonics in the 
vicinity of the resonant curve has began. Time t = 4370. 
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Figure 23: Decay of the monochromatic capillary wave. Zoom of the initial instability region. 
Time t = 1204r0. 
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7. Instability of the standing wave. 

May be the most practically important case of surface waves instability is the 
case of instability of the standing wave, when we have interaction of two waves 
a^^ and a_^^. In this case resonant curve is a circle with the center at zero wave 

number vector and of radius |fcol- It is clear that such a process is general for 
any isotropic dispersion. The theory for very similar instability in plasma was 
developed in [ 28,1 . 

7.1. Standing capillary wave. 

Simulation results for the standing capillary (/i = 0.1) wave are represented 
in Figures I25H29I Contour plots correspond to the level \a^^ = 10~^^, which is 
order of magnitude higher than background noise. The full picture of the k- 




Figure 24: Instability of the standing capillary wave. Initial conditions. Time t - 0. 

plane in the final moment of simulations is represented in Figure|30l Zoom of the 
most interesting region of the Figure [30] (initial instability region) is represented 
in Figure |3T] We observe isotropization of the wavefield, but in order to obtain a 
smooth spectrum we need to wait much longer. 
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Figure 25 : Instability of the standing capillary wave. Beginning of growth on the resonant curves. 
Some nonresonant absorption is noticeable. Time t = 14T0. 




Figure 26: Instability of the standing capillary wave. Growth on resonant curves continues. Time 

t = 57ro. 
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Figure 27: Instability of the standing capillary wave. Unstable harmonics are well developed. 
Timer = 2837^0. 




Figure 28: Instability of the standing capillary wave. One can notice formation of the forced 
harmonics as copies of initial circle shifted by pmkQ vectors. Time t - 509r0. 




Figure 29: Instability of the standing capillary wave, 
by secondary decay processes. Time t - lOlSrO. 



Isotropization started. 



The plane is filled 
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Figure 30: Instability of the standing capillary wave. FuU ^-plane. Time t - 2587r0. 
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Figure 31: Instability of the standing capillary wave. Zoom of the initial instability region. 
The inner circle corresponds to instability of the standing waves produced in the points of the 
maximum of the growth rate for decay of every individual initial waves. Compare position with 
the most developed harmonics inside the main circle in Figure|27] Time t - 2587r0. 
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7.2. Standing gravity wave. 

Simulation results for the standing gravity wave of steepness fi = OA are rep- 
resented in Figures l33]|37l Contour plots correspond to the level = 10"^^, 
which is order of magnitude higher than background noise. The full picture 




Figure 32: Instability of the standing gravity wave. Initial conditions. Time f = 0. 

of the fc-plane in the final moment of simulations is represented in Figure |38l 
Zoom of the most interesting region of the Figure ES] (initial instability region) 
is represented in Figure |39] Finally we observe almost complete isotropization 
of the wavefield, although we started from just two waves. Weak angle depen- 
dence resembles the cos(6) of coupling coefficient (|52l) . The observed process 
can be used for generation of isotropic wave field through initial generation of 
the standing wave, which in turn, through the discussed instability, will generate 
isotropic spectrum. This is quite nontrivial problem for direct wave generation. 

In our simulation we observed start of formation of the weakly turbulent 
spectrum tail (see Figure |40]) and formation of Kolmogorov-Zakharov weak tur- 
bulent spectrum of direct cascade [26] . 

In conclusion we have to note, that although the wave amplitudes should be 
high enough to make grid discreteness unimportant, at the same time they have 
to be low enough to satisfy weak nonlinearity conditions. 
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Figure 33: Instability of the standing gravity wave. Unstable harmonics on the resonant curves 
begin to grow. Time t = 1 1670. 
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Figure 34: Instability of the standing gravity wave. Unstable harmonics are well developed. 
Time t = 23270. 




Figure 35: Instability of the standing gravity wave. Formation of forced harmonics correspond- 
ing to initial circle shifter with ±ko vectors. Time t = 34870. 
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Figure 36: Instability of the standing gravity wave. Secondary processes reveal themselves. Time 
t = 463r0. 
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Figure 37: Instability of the standing gravity wave. Stochastization of the wavefield begins. Time 

t = 58oro. 
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Figure 39: Instability of the standing gravity wave. Zoom of the initial instability region. Time 
t = 306870. 
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Figure 40: Instability of the standing gravity wave. Beginning of the formation of the weakly 
turbulent spectrum tail. Solid line: angle averaged spectrum from simulation results; dashed line: 
theoretically predicted KZ-spectrum [26]; dotted line: the best fit with powerlike function. Time 
t = 3068rO. 
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8. Conclusion. 

We gave a complete theoretical descriptions of the three and four wave in- 
stabilities due to resonant interactions of waves. In order to simulate waves tur- 
bulence, this mechanisms have to work even on discrete grid, where exact reso- 
nance conditions are never fulfilled exactly. We demonstrated the possibility to 
achieve resonant interactions on discrete grid in numerical simulations. Simula- 
tion results are in good agreement with the theoretical predictions. 

Also we described in details algorithm for simulation of weakly nonlinear 
gravity-capillary surface waves. Numerical scheme, which was used in the code, 
conserves Hamiltonian of the system. Features of the algorithm are used to con- 
veniently control the adaptive time step. Described pseudo spectral method al- 
lowed us to simulate wave turbulence in numerous cases. 

We discussed and simulated instability of standing wave form both capillary 
and gravity waves. Numerical simulations show that instability of propagating 
waves leads to formation of anisotropic weakly-turbulent spectra, while insta- 
bility of standing waves leads to generation of almost isotropic spectra demon- 
strating tendency to formation of Kolmogorov-Zakharov tails. We conclude that 
numerical simulation of waves instability is a perfect tool for study of the wave 
turbulence theory. For experimental wave tanks this instability provides a very 
simple and robust approach which allows to produce isotropic wave field by ex- 
citation of just one standing wave. 

The paper as a whole can be used as a comprehensive guide for theoretical 
and computational approaches to simulation of weakly nonlinear waves on the 
surface of the fluid. 
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Appendix A. Discrete Hamiltonian variation. 

Let us derive variation of Hamiltonian (H"'^^ - H"). The k operator is self- 



adjoined 

J gkfd'r = j fkgd'r. 



(A.l) 



Let us perform variation in details for quadratic part of the Hamiltonian (fT3]) 
in the case of surface gravity waves. 



-'-'0 
1 r 

= - I {i^^^kf"^^ - ij/kilf"-^ + - (A"^'A")d^r+ 

I J (^""' - ri'^M^' + r]")d^r =^J [(i^"' - rfkr"' + 
.«,--,")]dV.f/(r-.'W-..")d^.= 

= i J A#(^""' + ^")dV + I J A77(77""^ + /7")dV = 

Here and further Ai/r = - i/^") and A?7 = (77"+' - 77"). 
Similar calculations give us all other variations. 
For short let us omit integral signs in varied expressions. 
Quadratic terms 

A|^J"#i/rdVj — > iA#(i/r"+^ +(A"); 



(1 



1 I ph/? 

2 J 



iA4(,-+,^). 



(A.2) 



(A.3) 



Cubic terms 
A 



-|A(A(V,(77"+i + ?7")V((A"+i + (A")) 

+|A77(|V^«+i|' + |V^"p); 
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+ 



(A.4) 



Quartic terms 



X k(rf^^kif/"+^ +ri"ktf/")] 
xk(rf'*^kil/"*^ +7/4iA")]; 



+ 



{Til/){kil/)ifd^r 



+ 



+|A#[((77"+i)2 + (?7")2) X 
X V2(i/r«+l + + 
+|A?7(?7"+' + ?7")X 



Appendix B. Matrix elements 



We repeat formulae from [|23ll . 

4;rV2 (\5;tA^,A;tJ 



BhAv, B 



1/4 



1 



AkBki Bk2 



1/4 



L^'\kuk2) 



= 1^1, h=g + o-k^. 
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y(2,2) _ 2 
k,ki ,k2,k 



a? 

k 



-2 J I<\k2Xki)blbj^^6(k + k2-k2)dhdk2 
+ f r^^\t Jcu li2)bl b* 6(Jc + lci+ 4)d^id4, 

J h k2 

a^^^ = J B(k, kiM, h)blb^b^^5{k -kx-h- 4)d^id4d4 



2oJk- oJh - (^k2 ' 



T<\iA,k2) = -\ 



2ojt- coki - (i>k2 ' 

- Y^\t 4, t - 4)r^'^(4, ^1, 4 - ^i) 

- r^^^^ +ki,t ti)T^^\t2 + 4, k, ici) 

+ r(2\-l -4,4 4)r^'^(-4 - 4, 4, 4). 
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7^1234 - 2 (-^1234 + 7^2134)> (B-H) 
1 1 

^^^^^ ~ ~ lenHhkjk^hy/'^ 

-l2k\k2kj,kA - 2(a»i + 0*2) [a»3a»4((fci • ^2) - ^1^2) 
1 

+ oJiOJiiih • ka.) - hk^)] — 

r 



- 2((x»i - a>4) [a)20J^{{ki ■ ki) + ^1^4) + a)y(jJi^{{k2 ■ ^3) + ^2^3)] — 



1 



+ [(^1 • ^2) + ^1^2] [(^3 • h) + ^3^4] 

+ [-(^1 ■ ^3) + ^1^3] [-(4 ■ 4) + ^2^4] (B.12) 

+ [-(^1 ■ h) + hkAW-ih ■ h) + ^2^3] 

2 [(^1 ■ ^2) - ^i^2][-(4 ■ 4) - ^3^4] 
+ 4(a>i + C02) -, ■ 



2 [(^1 ■ ^3) + kihMih ■ ki) + ^2^4] 



2 [(^1 • h) + ^1^4][(^2 • ^3) + ^2^3] 



Appendix C. On the free surface hydrodynamic model 

In details this question was considered in ^sT] but here we shall follow orig- 
inal consideration which was done by A. I. Dyachenko in 1995 (result was men- 
tioned in [32]). It is shown that adequate free surface hydrodynamic model can 
be obtained when taking into account at least fourth order term in the Hamilto- 
nian written for variables rj and if/. 

1. Let's consider the fourth order Hamiltonian in one dimensional case. It 
has a form: 



Hilj, + iJ/H{r]H{r]Hifi,YY}dx+ 
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H = 

2 J ' ' 
^ 2 



+ ^ J" {gn^ + o-ril}dx 

Here H is the Hilbert operator. The equations of motion 

_6H _ 6H 



and 



— = - //(A. - H(r]Hi^,), - (77(Ax).- 



(5(A 

SH 1 

— = - -((^lA.)' - ^l) - ^x.riHip, - Hip, X H{TjHi/f,),+ 
orj 2 

Now let's consider the small scale perturbation Sif/ and 6ri on the large scale 
background i^o and 770 

Perturbations 6^ and 6^ acquire the form: 
SH 

6— = kdip - (ktPo)k6r] - iVi//o)^dr] 
01]/ 

6H , , 

07] 

(Note: terms proportional to 770 and V^6t] vanish except crV^S?].) 
From (IC.ll) one can get the following dispersion relation: 

iaj,-kvof=g\k\ + (T\k\' (C.l) 

When taking into account only cubic Hamiltonian the {kif/ofkdri term in (|C.1I) is 
absent and dispersion relation is 'bad' (instability at large k): 

{Ojk - kvof = g\k\ + (T\k\' - {klfjofk^ 
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